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ABSTRACT 

We assess the impact of inhomogeneous reionization on detection of HII regions sur- 
rounding luminous high redshift quasars using planned low frequency radio telescopes. 
Our approach is to implement a semi-numerical scheme to calculate the 3-dimensional 
structure of ionized regions surrounding a massive halo at high redshift, including the 
ionizing influence of a luminous quasar. As part of our analysis we briefly contrast our 
scheme with published semi-numerical models. We calculate mock 21 cm spectra along 
the line of sight towards high redshift quasars, and estimate the ability of the planned 
Murchison Widefield Array to detect the presence of HII regions. The signal-to-noise 
for detection will drop as the characteristic bubble size grows during reionization be- 
cause the quasar's influence becomes less prominent. However, quasars will imprint a 
detectable signature on observed 21 cm spectra that is distinct from a region of typical 
IGM. At epochs where the mean hydrogen neutral fraction is » 30% or greater we 
find that neutral gas in the IGM surrounding a single quasar will be detectable (at a 
significance of 5cr) within 100 hour integrations in more than 50% of cases. 1000 hour 
integrations will be required to detect a smaller neutral fraction of 15% in more than 
50% of cases. A highly significant detection will be possible in only 100 hours for a 
stack of 10 smaller 3 proper Mpc HII regions. The accurate measurement of the global 
average neutral fraction ((xhi)) will be limited by systematic fluctuations between 
lines of sight for single HII regions. We estimate the accuracy with which the global 
neutral fraction could be measured from a single HII region to be 50%, 30% and 20% 
for (xhi) ~ 0.15, 0.3 and 0.5 respectively. 
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1 INTRODUCTION 

The reionization of cosmic hydrogen is commonly believed 
to have been due to UV photons produced by the first stars 
and quasars and wa s an important milesto ne in the history 
of the Universe (e.g. lBarkana fc Loebll200ir h The recent dis- 
covery of very distant quasars has allowed detailed absorp- 
tion studies of the state of the high redshift intergalactic 
medium (IGM) at a ti me when the Un i verse was less tha n 
a billion years old (e.g. iFan et alj|2005 IWhite et aLll2003h . 
Several studies have used the evolution of the ionizing back- 
ground inferred from these spectra to argue that the reion- 
ization of cosmic hydrogen was completed just beyond z = 6 
jFan et all 120061 : iGnedin fc Fanl 120061 ; IWhite et al.l l2003h . 
However, other authors have claimed that the evidence for 
this rapid change is the result of an inc orrect choice for th e 
distribution of overdensities in the IGM l|Becker et alj20"07h . 

Several of the most distant quasars exhibit complete 
Gunn-Peterson troughs, the red edges of which do not ex- 



tend as far as the redshifted Lycv wavelength. The sim- 
plest interpretation of this feature is the presence of an 
HII r egion in a partially neutral IGM (e.g. ICen fc Haimanl 
2000). This idea has been used to argue that there is a sig- 
nificant fraction of neutral hy d rogen in the IGM beyond 
z sa 6 JWvithe fc Loebl l2004al ; iMesinger fc Haimanl 12004 
IWvithe et aTT l2005). with a lower limit on the hydrogen 
neutral fraction that is as many as two orders of mag- 
nitude larger than constraints provided by direct absorp- 
tion studies. On the other hand, more detailed numerical 
modeling h as sho wn that this interpretatio n is uncertain 
|Lidz et all 2007a; iBolton fc Haehneld l2007h , and that the 
observed spectra could either be due to an HII region or a 
classical proximity zone. 

The difficulties in interpreting spectra of high redshift 
quasars arise because the large optical depth for Lya ab- 
sorption in a neutral IGM limits its effectiveness in probing 
reionization at redshifts beyond z m 6. A better probe of the 
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reionization history involves the redshifted 21 cm emission 
from neutral hydrogen in the IGM. Several probes of the 
reionization epoch in redshifted 21 cm emission have been 
suggested. These include observation of the emission as a 
function of redshift averaged over a large area of sky. This 
observation would provide a direct probe of the evolution in 
the neutral fraction of the IGM, and is referred to as the 
global step i|Shaver et al.lll999l ; iGnedin fc Shaved I2004I ). A 
more powerful probe will be provided by observation of the 
power spectrum of fluctuations together with its evolution 
with redshift. This observation would trace the evolution 
of neutral gas with reds hift as well as the topology of the 
reionization process (e.g. Tozzi et a l. 2000; Furlanetto et al.l 
2004al; lLoeb fc Zaldarriagal |2004 IWvithe fc Morales! 120071; 



Iliev et all 20061 ). Finally, observation of individual HII re- 



gions will probe quasar physics as wel l as the evolution o f 
the neutral gas (IWvithe fc Loeb|[2004bl ; iKohler et ai]|2005h . 
iKohler et~a l, (2005) have generated synthetic spectra using 
cosmological simulations and conclude that quasar HII re- 
gions will provide the most prominent individu al cosmolog- 
ical signals. Recent work bv lDatta etail (|2007t ) has focused 
on the detection of ionized bubbles in redshifted 21 cm maps 
using a visibility based formalism. Their results suggest it 
may be possible to blindly detect spherical HII regions of 
radius J> 22 comoving Mpc during the epoch of reionization. 

Various experiments are planned to measure 21 cm 
emission from the pre-reioniz ation IGM, includin g the Low 
Frequency ArrajQ (LOFAR: IValdes et ail l2006h and the 
Murchison Widefield ArrajQ (MWA). The MWA is being 
constructed on a very radio quiet site in Western Australia 
and will consist of 512 tiles, each containing 16 cross dipoles, 
which operate over the range 80-300 MHz. The tiles will be 
spread to cover baselines of up to 1.5 km. In collecting area 
the MWA would be comparable to the Very Large Array 
(VLA). However, the correlator/receiver system will be able 
to handle a much larger number of channels than are avail- 
able on the low frequency VLA. In this paper we focus on the 
observation of individual HII regions around known quasars 
with emphasis on the capabilities of the MWA. In particular 
we address the influence of the percolating IGM on the red- 
shifted 21 cm signal from the neutral hydrogen surrounding 
a quasar HII region. Our analysis extends ideas presented by 
IWvithe et al.l (|2005) by including calculations of a realistic 
structure for ionization of the IGM, as well as estimates of 
the achievable signal-to-noise. 

We discuss a simple detection strategy consisting of a 
spectral fit to the 21cm signal along the line of sight to- 
wards a known high redshift quasar. The following back- 
et f -the- envelope calculation gives a comparitive gauge of the 
statistical significance required by this strategy relative to a 
blind search for individual HII regions, embedded in a par- 
tially ionized IGM. 

For a blind search, the number of independent samples 
within a single MWA field is given by 

Vmwa Vmwa 



N s = 



ttRI x 27? q ' 



where Vmwa is the comoving volume of the MWA field and 
V q (z, 7? q , Rh) is the intersection of the synthesized beam (of 



1 http://www.lofar.org/ 

2 http:/ /www. haystack. mit.edu/ast/arrays/mwa/ 



comoving radius Rh) with the HII region (of comoving radius 
7?q). For a 4.5 proper Mpc region at z = 6.5 using a synthe- 
sized beam of width (9b = 3.2' we obtain V q ~ 10 4 comoving 
Mpc 3 . Together with Vmwa ~ tt4700 2 x 260 ~ 10 10 comoving 
Mpc 3 (|Bowman et al.ll2006h . we have iV s ~ 10 6 . 

Given a signal-to-noise level of na for detection of the 
mean IGM, the number of false positive detections per MWA 
field is 



Afaiso ~ iV s erfc 



( — 

\V2a 



(2) 



where erfc(x) is the complimentary error function. For exam- 
ple, if the mean IGM can be detected with a signal-to-noise 
level of 2a, we expect ~ 10 4 false detections per field, but 
only one true detection of a quasar HII region. The chance 
of obtaining only one false detection requires n ~ 5. There- 
fore, in order to reduce the chance of a false detection and 
have confidence in the reality of the HII region we require 
Afaiae <C 1 or n 2> 5. This level of significance will probably 
not be possible with the expected sensitivity of the MWA 
which will make a blind search impractical. 

On the other hand, by making an observation around 
a quasar, we are able to take a Bayesian approach by uti- 
lizing the likelihood of detection given that an HII region 
is known a priori to be present. In this case a detection at 
a given significance level is not diluted by the possibility of 
false detections. We stress that even when the impact of the 
quasar is no longer dominant, the method of detecting the 
21 cm signal a posteriori along the line of sight to a known 
quasar will still succeed, because neutral gas is known a pri- 
ori to be absent in a large volume surrounding the quasar. 

In addition to its physical configuration, the sensitivity 
of a radio telescope is limited due to difficulties in achiev- 
ing a perfect calibration and to problems introduced by 
foregrounds and the ionosphere. There has been significant 
discussion in the literature regarding foregrounds. While 
fluctuations due to foregrounds are orde rs of magnitude 
larger than the reioni zation signal (e.g. |Pi Matteo et al.l 
2002: 1 Oh fc Ma ck 2003), foreground spectra are anticipated 
to be smooth. Since the reionization signature includes fluc- 
tuations in frequency as well as ang le, they can be removed 
through con t inuum subtraction (e.g. lGnedin fc Shaverll2004l : 
IWang et all l200d ). or using the differences in symme- 
try from power spectra analysis (M orales fc Hewitt] 120041 : 
IZaldarriaga et akl 120041 '). In this paper we concentrate on 
theoretical limitations of the MWA and the impact of the 
inhomogeneous reionization of the IGM. The practical dif- 
ficulties introduced by calibration and foreground removal 
etc. need to be addressed by full simulations of interfer- 
ometric observations, and ultimately with early data sets 
from new facilities. 

We begin by describing the inclusion of stellar reioniza- 
tion and a luminous quasar within a semi-numerical scheme 
to produce realistic fluctuations in the ionized gas distribu- 
tion (§[2|. We then use this model to describe the evolution 
of the ionization structure of the IGM surrounding a mas- 
sive quasar host (both with and without the presence of a 
quasar HII region) in §[3l before discussing the detectability 
of these regions in §[4] We present our conclusions in §[5] 
Throughout this paper we adopt t he set of cosmologic al pa- 
rameters determined by WMAP (jSpergel et al.| [2007) for a 
flat ACDM universe. 
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2 SEMI-NUMERICAL MODEL FOR THE 

GROWTH OF QUASAR HII REGIONS IN A 
PERCOLATING IGM 

In this section we describe our model for the reionization of 
a 3-dimensional volume of the IGM using a semi-numerical 
scheme tha t is sim ilar in spirit to the models descr ibed by 
IZahn et al, I (120071') andlMesinger fc Furlanettol (|2007l ). based 



Furlanetto et all (|2004bF ~ We begin with a brief de- 



upon 

scription of the generation of the Gaussian initial conditions 
in the density field surrounding a massive quasar host halo 
at z ~ 7. We then describe a semi-analytic model for the 
density dependent reionization of the IGM, including the in- 
fluence of a nearby quasar, before combining the numerical 
realization of the density field with our semi-analytic calcu- 
lations to generate 3-dimensional realizations of the ioniza- 
tion state of the IGM. 



2.1 Construction of density field surrounding a 
massive halo 

We simulate the linear density field <5(k) inside a periodic, 
comoving, cubic region of volume V = L 3 , by calculating 
the density contrast in Fourier space S(k) corresponding to 
a ACDM power spectrum at a specified redshift. The power 
spectrum used has been normalized using as and is given by 
P(k,z) oc k"T 2 (k)D 2 (z) where T(k) is an an alytic approxi- 
matio n of the present day transfer function l|Bardeen et al.l 
ll986T). n = 0.96 is the primordial power spectral index 
I Spergel et ail 12007ft and D{z) is the linear growth factor 
between redshift z and the present. 

Numerically, this procedure is carried out by laying 
down a cubic grid of size N 3 and approximating the con- 
tinuous Fourier relationship between 5(x) and <5(k) with its 
discrete form, 



N 3 



(3) 



Each Fourier density mode is calculated by sampling the 
power spectrum at each of the finite set of wavenumbers k = 
(I, m, n)k , where I, m, n € {-N/2 + 1, -N/2 + 2,..., N/2} 
and ko = 2n/L. In order to realize the Gaussianity of 
<5(x) we draw both the real and complex components of 
<5(k) from a zero-mean normal distribution with variance 
a 2 = P(k, z)V/2. Furthermore, we enforce Hermiticity upon 
<5(k) to ensure 5(x) is real. All 7 vertices corresponding to 
the positive frequency Nyquist modes must be real, and in 
order to produce a zero-mean density field the DC mode 
<5(k = 0) is set to zero. The 3-dimensional inverse Fourier 
transform of S(k) gives an Eul erian realiz ation of <5(x) fil- 
tered on a scale A = L/N (e.g. ISirkol [20051 '). 

To study the impact of a high redshift quasar on the 
ionization state of its surrounding IGM we identify the 
approximate location of a collapsed halo of mass M by 
smoothing the linear density field using a real-space spher- 
ical top-hat filter of radius _Rf = (3M/47rpo) 1//3 , where 
po = 3f2 m //o/87rG is the present mean matter density. This 
is equivalent to multiplying the Fourier coefficients of S(k) by 
the normalized kernel W(k) = 3ji(kR{)/kR{, where ji(x) is 
the first-order spherical Bessel function. Density peaks are 
considered to have collapsed if <5(x, z) > S c , where 8 C w 1.686 
is the critical linear overdensity of a spherical top-hat den- 



sity perturbation. Note that the field is linearly evolved to z 
before the <5 C filter is applied. We fiducially set M = 10 12 Mq 
as the minimum value for the masses of quasar host halos 
at z 6.5 with proximity zones of R w 30 comoving Mpc 
|Fan et alj|2006h . 



2.2 Density dependent model of global 
reionization 

Our semi-analytic model to compute the relation between 
the local dark matter overdensity and the re ionization of the 
IGM is bas ed on the model describe d by IWvithe fc Loebl 
(|2007h and IWvithe fc Morales! (|2007f ). Here we summarize 
the main features of the model and refer the reader to those 
papers for more details. 

Any model for the reionization of the IGM must de- 
scribe the relation between the emission of ionizing photons 
by stars in galaxies and the ionization state of the inter- 
galactic gas. This relation is non-trivial as it depends on 
various internal parameters (which may vary with galaxy 
mass), such as the fraction of gas within galaxies that is 
converted into stars and accreting black holes, the spectrum 
of the ionizing radiation, and the escape fraction of ionizing 
photons from the surrounding interstellar medium as well as 
the ga lactic halo and its immediate infall region [see lLoebl 
(2006) for a review]. The relation also depends on inter- 
galactic physics. In overdense regions of the IGM, galaxies 
will be overabundant because small-scale fluctuations need 
to be of lower amplitude t o form a galaxy whe n embedded in 
a larger scale overdensity (|Mo fc Whitdfl996T ) . On the other 
hand, the increase in the recombination rate in overdense 
regions counteracts this galaxy bias. The process of reioniza- 
tion also contains several layers of feedback. Radiative feed- 



galaxy formation (lEfstathioulll992 


Thoul & Weinberg! 


iQuinn et al.l Il996: Diikstra et al. 


2004). This delays 



the completion of reionization by lowering the local star for- 
mation rate, but the effect is counteracted in overdense re- 
gions by the biased formation of massive galaxies. 

The evolution of the ionization fraction by mass Qs,r of 
a particular region of scale R with overdensity 8 (at observed 
redshift z bs) may be written as 



dQs,R 
dt 



N io 



0.76 



!s,r- 



dF col {8,R,z,Mi 

O 



dt 

dF co i(5,J?,2,M min 



D(z) 



dt 



(1+*) 



16,R, 



(4) 



D(Zobs) 

where N[ on is the number of photons entering the IGM per 
baryon in galaxies, cub is the case-B recombination coeffi- 
cient, C is the clumping factor (which we assume, for sim- 
plicity, to be a constant value of 2) and D(z) is the growth 
factor between redshift z and the present. The production 
rate of ionizing photons in neutral regions is assumed to be 
proportional to the collapsed fraction F co i of mass in ha- 
los above the minimum threshold mass for star formation 
(M m i n ), while in ionized regions the minimum halo mass 
is limited by the Jeans mass in an ionized IGM (Mi on ). 
We assume M m i n to correspond to a virial temperature 
of 10 4 K, representing the hydrogen cooling threshold, and 
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Miod to correspond to a virial temperature of 10 s K, repre- 
senting the mass below which inf all is suppressed from an 
ionized IGM (|Diikstra et al.l l2004). In a region of comoving 
radius R and mean overdensity S(z) — SD(z) / D(z hs) (spec- 
ified at redshift z instead of the usual z = 0), we use the 
extended Press-Sc hechter l|Press fc Schechterl Il974h model 
(Bon d~t al .11 199 il l to calculate the relevant collapsed frac- 
tion 

F col (5, R, z) = erfc [ ^ZJM -_ ) , (5) 

where erfc(:r) is the complimentary error function, cr 2 (R) is 
the variance of the density field smoothed on a scale R, and 
o-g al is the variance of the density field smoothed on a scale 
Rgai, corresponding to a mass scale of Af m j n or Afj on (both 
evaluated at redshift z rather than at z = 0). 

Equation @ may be integrated as a function of 8. At 
a specified redshift this yields the filling fraction of ionized 
regions within the IGM on various scales R as a function 
of overdensity. We find that this model predicts the sum of 
astrophysical effects to be dominated by galaxy bias, and 
that as a result overdense regions are reionized first. This 
leads to the growth of HII regions via a phase of percolation 
during which individual HII regions overlap around clustered 
sources in overdense regions of the universe. We may also 
calculate the corresponding 21 cm brightness temperature 
contrast 

T(6,R) = 22 mK (1-Qs,r) H^"J i 1 + *) ■ (6) 

In this paper we ignore the enhancement of the bright- 
ness temperature fluctuations due to peculiar velocities in 
overd ense regions (|Bharadwai fc A li 2005; Barka na fc Loebl 
2005). Peculiar velocities were included in the semi- 
numerical model of iMesinger fe Furlanettol ([2007) who 
found their effect to be small on scales ~ 10 comoving Mpc. 

2.3 The ionization field 

Having determined the value of the ionized fraction Qs,r 
as a function of overdensity 5 and smoothing scale R, we 
may now construct the ioniza tion field. We employ a s imilar 
filtering algorithm to that of IMesinger fc Furlanettol (120071 ) 
to determine the ionization state at each grid position. This 
is done by repeatedly filtering the linear density field on 
scales in the range from L to L/N at logarithmic intervals 
of width ARt/Rt = 0.1. For all filter scales, the ionization 
state of each grid position is determined using Qs,r and 
deemed to be fully ionized if Qs,r ^ 1. All voxels within a 
sphere of radius R centered on these positions are flagged 
and assigned Qs,r = 1, while the remaining non-ionized 
voxels are assigned an ionized fraction of Qs,R !mill , where 
^f,min = L/N corresponds to the smallest smoothing scale. 
A voxel forms part of an HII region if Qs,r > 1 on any scale 
R. 

2.4 Inclusion of quasars in the semi-numerical 
scheme 

In a region of radius R containing a quasar, the cumulative 
number of ionizations per baryon will be larger than pre- 



dicted from equation Q. We include quasars in our scheme 
by first computing the fraction of the IGM within a region 
of radius R centered on the position x that has already been 
reionized by stars [using equation (J2J] , we then add an addi- 
tional ionization fraction equal to the quasar's contribution 

where R q is the radius of the HII region centered on x q that 
would have been generated by the quasar alone in a fully 
neutral IGM. In contrast to stellar ionization the quasar 
contribution comes from a point source, and so the contri- 
bution to Q originates from a single voxel only (rather than 
all voxels within |x — x q |). Following this addition we filter 
the ionization field as described in § 2.3. 



3 THE EVOLUTION OF HII REGIONS IN A 
PERCOLATING IGM 

In this section we describe one example of the evolution 
of the IGM during the reionization epoch. We first show 
a model for the evolution of ionized regions due to stellar 
ionization alone in §3.1, before adding the presence of a 
luminous quasar in § 3.3. 

3.1 Fiducial model for reionization 

Throughout this paper we consider a model in which 
the mean IGM is r eionized at z = 6 (|Fan et al ] 120061 ; 
iGnedin fc Far] 120061 ; I White et all 120031 '). In this model we 
assume that star formation proceeds in halos above the hy- 
drogen cooling threshold in neutral regions of IGM. In ion- 
ized regions of the IGM star formation is assumed to be 
suppressed by radiative feedback (see §[2]). We compute the 
3-dimensional structure for this model using the prescrip- 
tion outlined in §[2j In the left-hand panel of Figure [T] we 
present snapshots of the ionization fraction Q(x) for slices 
through a realization of this model when the mean, volume- 
weighted neutral fraction of the IGM is « 0.5, 0.3, 0.15 and 
0.1 corresponding to redshifts z — 7.5, 7.0, 6.7 and 6.5 re- 
spectively. The simulation presented corresponds to a linear 
density field of resolution 256 3 , with a comoving side length 
of 219 hT 1 Mpc. Each slice is 0.86 h" 1 comoving Mpc deep. 
Since the box is centered on a massive halo the simulation 
shows a large HII region growing in the center of the field. 
As reionization progresses we see this central region merge 
with adjacent HII regions. This proce ss of overlap is seen in 
detailed numerical simulations (e.g. iMcQuinn et all I2007T) 
and leads to an irregular structure of neutral gas filaments 
late in the reionization era. The overlap of some HII re- 
gions, combined with the fact that our model does not allow 
for the ionizing flux from one source to contribute to the 
ionization due to another results i n a violation of photon 
conservation (e.g. IZahn et alj|2007h . Although this effect is 
expected to be small during the early stages of reioniza- 
tion when there are fewer HII region mergers, it does lead 
to more severe discrepancies at later times when the global 
ionization fraction approaches unity. Figure [2] illustrates this 
effect by plotting the global neutral fraction as a function 
of redshift from our simulations. The dashed line shows the 
evolution of mean neutral fraction by mass as expected from 
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Figure 1. Ionization maps of the IGM showing the evolution of global ionization fraction, Q, with and without the ionizing contribution 
of a quasar. The slices shown correspond to epochs when the global, volume-weighted neutral fraction is (^hi) ~ 0.46, 0.28, 0.16 and 
0.08 (top to bottom) corresponding to z = 7.5, 7.0, 6.7 and 6.5 respectively. Each slice has a comoving side length of 219 ft -1 Mpc and 
is 0.86 ft. -1 comoving Mpc deep. The ionization maps were computed within boxes of dimension 256 3 . Left: Evolution with no quasar. 
Middle left: Corresponding ionization maps with an embedded R q ~ 34 comoving Mpc quasar HII region. Dashed circles represent the 
HII region that would have been generated by the quasar alone in a fully neutral IGM. Middle right: Corresponding synthetic maps with 
$b ^ 3.2' (corresponding to central peak to first null A9 w 5.5') including quasar HII region using the MWA and 1000 hr integrations. 
Right: Corresponding synthetic maps with 6*t, pa 5.8' (corresponding to central peak to first null A8 w 10.0') including quasar HII region 
using the MWA and 1000 hr integrations. 



the analytic model which provides the input density depen- 
dent ionization fraction. We find that our semi-numerical 
scheme yields results that conserve photons at the ~ 15% 
level for neutral fractions (xjji) <; 0.2. 

3.2 Comparison with other semi-numeric models 

We note that our scheme correctly assigns the ionization 
fraction on scales that are not resolved by the simulation 
cube. This is not necessarily the case in similar models that 
approximate the ionization field as a two-phase medium (e.g. 



iMesinger fc Furlanetto|[2007l ) , although a two-phase approx- 
imation may be justified depending upon the resolution of 
the simulation. Implementing our scheme with increased res- 
olution would improve the accuracy and dynamical range 
of our simulations. However, we maintain a high resolu- 
tion relative to the synthesized beam. The assignment of 
regions as either fully neutral or fully ionized may also lead 
to an overestimate of the brightness temperature contrast 
across ionization fronts which are observed at finite resolu- 
tion. Our implementation achieves approximately the cor- 
rect relation between neutral fraction and epoch over most 
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Figure 2. Global, mass- weighted neutral fraction as a function of 
redshift. The dashed line shows the evolution of neutral fraction as 
expected from the density dependent model described in § 2.2 (for 
an overlap at z 6). The data show the mass- weighted globally 
averaged neutral fraction from our simulations. Hollow squares 
indicate the approximate epoch for which simulations were used 
to obtain the mock line of sight spectra shown in Figure 111 

of the reionization history and produces soft transitions at 
the edges of large HII regions by smoothing over HII re- 
gions not resolved by the sim ulations. However, as noted by 
iMesinger fc Furlanettol ([20071 ) a two-phase approximation is 
appropriate for the purposes of describing the morphology 
and evolution of HII regions. We also point out that there 
are several features not incorporated in our model. These 
include the galaxy distribution (which will enhance fluctu- 
ations on small scales due to Poisson noise) and the veloc- 
ity field of the gas (see IMesinger fc Furlanettoll2007l ). How- 
ever these additional considerations are not important on 
scales larger than 10 comoving Mpc, and so do not influence 
our results regarding HII regions surrounding high redshift 
quasars. 



has a comoving side length of 219 /i" 1 Mpc and is 0.86 
comoving Mpc deep. Note that while the evolution of the 
ionization due to stars has been followed self-consistently in 
these simulations, we have assumed the quasar to contribute 
the same proper volume of ionization at each redshift. The 
apparent evolution in the size of the HII region surrounding 
this quasar is therefore due to the additional contributions 
from nearby galaxies to the HII region (and not due to the 
evolution of IGM density since the slices are shown in co- 
moving units). 

The quasar HII region appears spherical and distinct 
early in the reionization era when the global neutral frac- 
tion is still large, and the typical size of galaxy driven HII 
regions is small (i.e. z > 8 in this model). However, the 
shape of the HII region becomes highly distorted once the 
typical galaxy (or galaxy cluster) generat ed HII region be- 
comes large late in the reionization epoch dFurlanetto et all 
l2004ar ). Indeed, once the size of the typical HII region ex- 
ceeds that generated by the quasar, the quasar HII region 
ceases to be an identifiable feature. We note that a band- 
width of Av, centered on z, corresponds to the redshift in- 
terval Az = Ai/(1 + z) 2 /1420 MHz, implying that a 20 MHz 
band (the approximate frequency range of our synthetic 
spectra) has a redshift extent of Az m 0.8. We expect signif- 
icant evolution of the 21 cm signal over this redshift interval 
very close to overlap. However, given the other uncertainties 
in our modeling we ignore such evolution over this band- 
width in our simulations. 



4 DETECTABILITY OF HII REGIONS 

In the remainder of this paper we discuss estimates for the 
detectability of the predicted signature of quasar HII regions 
in 21 cm emission. In particular, we discuss the sensitivity 
of the MWA to the HII regions surrounding known high 
redshift quasars in a percolating IGM. 



3.3 The influence of a quasar in the fiducial 
model for reionization 

We next add quasar luminosity to the simulations shown 
in the left-hand panels of Figure [T] using the prescription 
outlined in § 12.41 The value of _R q is subject to large un- 
certainties, including the quasar duty cycle, lifetime and lu- 
minosity. Assuming a line of sight ionizing photon emission 
rate of A r 7 /47r photons per second per steradian, a uniform 
IGM and neglecting recombinations, the line of sight extent 
of a quasar's HII region (from the quasar to the front of the 
region) observed at time t is (|White et al.ll2003T ) 



R q ~ 4 pMpc x m 



1/3 



10 57 10 7 yr J 



1/3 



l + Z 

7.5 



(8) 



where xm is the neutral fraction and t asc — t — _R(t ag e)/c 
is the age of the quasar corresponding to the time when 
photons that reach R q at time t were emitted. 

We show examples of R q « 34 comoving Mpc, consis- 
tent with the lower limit a round the luminous SDSS quasars 
at z > 6 (|Fan et alj 2006). The second column of Figure [1] 
shows the same realization of the IGM presented in the first 
column, but now includes the influence of a luminous quasar 
hosted by the central massive halo. As before, each slice 



4.1 Sensitivity to the 21 cm signal 

Radio interferometers measure a frequency-dependent, com- 
plex visibility V(\J,v) (Jy) for each frequency channel and 
baseline U in their configuration. The measured visibility 
is, in general, a linear combination of signal, foreground and 
system noise, 



V(U,u) = Vs(U,i/) + Vf(U,v) + Vn(U,u) 



(9) 



where Vs is the signal, Vf the contribution due to foreground 
sources and Vs the system noise. The signal of interest in 
this work is redshifted 21 cm emission. In order to discuss 
the response of a phased array to the brightness temperature 
contrast of the 21 cm emission from the IGM we must i) 
quantify the instrumental noise-induced uncertainty in each 
visibility AVn and ii) address foreground subtraction. The 
rms noise fluctuation per visibility per frequency channel is 
given by 

2/CB T$y S 



AVk 



(10) 



where T sys is the system temperature (K), A a g the effective 
area of one antenna (m 2 ) , t\j is the integration time for that 
visibility (s), Av is the frequency bin width (Hz) and fee 
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is the Boltzmann constant. Instrumental n oise is uncorre- 
cted in the frequency domain. As shown bv lMcQuinn et all 
(2006), the average integration time tu that an array ob- 
serves the visibility U is 



kj. [Mpc' 1 



0.2 



tu 



■n(U), 



(11) 



where tint is the total integration time and n(U) is the 
number density of baselines that can observe the visibil- 
ity U. The transverse wavenumber k± is given by k± = 
27rjU|/r cm (,z), where r em (z) is the proper distance to the 
point of emission. 

The MWA will consist of 512 tiles, each with 16 
cross-dipoles with 1.07 m spacing and an effective area 
Aeff ~ 16(A 2 /4) for A <; 2.1m. The system temperature 
at v < 200 MHz is sky domi nated and has a value T Bys ~ 
250 [(1 + z)/7] 2 ' 6 K. Following [Bowman et aD (|2006T ). we as- 
sume a smooth antenna density profile p ant oc r~ 2 within 
a 750 m radius with a core density of one tile per 18 m 2 . 
Due to the large system temperature, weak signal and poor 
sampling of the largest baselines, naturally weighted visi- 
bility data produce, at best, low signal-to-noise images and 
spectra. In addition, finite uv-coverage, determined by the 
maximum baseline, truncates the simulated visibility data 
and hence we are unable to use a Gaussian tapering func- 
tion to suppress small scale sidelobes. In order to maintain 
the highest sensitivity for spectral fitting (see §4.3), we in- 
stead truncate the naturally weighted visibility data using 
a filled circular aperture of radius a with a corresponding 
beamwidtrjf] (full width at half maximum) Oh = 0.705?7max, 
where (7 ma x = a/A. Figure shows the baseline density pro- 
file for the MWA configuration described in §4.1 at v « 
165 MHz, as well as the rotationally invariant system noise 
per visibility showing truncated visibilities due to variation 
of synthesized beamwidth. The third and fourth columns of 
Figure [1] show the resulting synthetic 21 cm brightness tem- 
perature maps of the quasar generated HII regions shown in 
the second column for beamwidths #b w 3.2' and 5.8'. For 
the purposes of imaging, uniform weighting of the visibilities 
would remove the bias from the greater sampling of short 
baselines (i.e. low spatial frequencies) and therefore lessen 
the masking of small scale features in the image. 

We simulate the thermal noise in a 3-dimensional 
visibility-frequency cube using equations (|10|) and pip . This 
procedure allows us to account for a particular array con- 
figuration and observation strategy. We then perform a 2- 
dimensional inverse Fourier transform in the wu-plane for 
each binned frequency in the bandwidth, which gives a re- 
alization of the system noise in the image cube (i.e. sky 
coordinates). Using the linearity of both the Fourier trans- 
form and equation ((SJ, the observed specific intensity can 
then be found either by first adding visibility terms and 



3 Note that the 'beamwidth' referred to by many authors is 
A9 = 1.22A/D ~ A/D, which is the highest angular resolution 
meeting the Rayleigh criterion. This corresponds to the central 
peak to first null separation in the resulting beam pattern of 
a circular aperture. Other conventions include: full beam width 
between first nulls $fwfn = 2.44A/Z5, full beam width at half 
power Ofwhp = 1.02A/D and full beam width at half maximum 
$FWHM = 0.705A/D. We have adopted the last of these conven- 
tions in this paper. 
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Figure 3. Top: Baseline density profile for the MWA configu- 
ration described in §4.1 at v ss 165 MHz. Bottom: Rotationally 
invariant system noise per visibility showing truncated visibili- 
ties (shaded regions) due to variation of synthesized beamwidth 
(#b = Sfwhm = 0.705t/mit f° r a circular aperture). 



then transforming to image space, or by adding the real- 
ization of each visibility component in image space. Finally, 
we apply the Rayleigh- Jeans approximation to express the 
observed flux in terms of brightness temperature T using 
dByjdT = 2/cb/A 2 , where £>„ is the specific intensity. 

The resulting rms noise in an image constructed in the 
manner described, for a beamwidth 8^ and frequency bin 
Av has the form 



AT b w T(0b) 



1 



Av 



100 kHz 100 hr 



-1/2 



(12) 



We find that the prefactor takes values of T{9 D ) ~ 11 mK 
for 8h ~ 3.2 ' (corresponding to A6 = 5.5 ' central peak to 
first null) and T(6 h ) w 65 mK for 9 b = 2.9 ' (AO = 5.1') 
For c omparison, the radiometer equation (e.g. IWvithe et al.l 
|2005( ) estimates the noise within a synthesized beam of width 
9h to have corresponding values of T(6h) ~ 20 and 22 mK 
for 9h = 3.2 ' and 2.9 ' respectively. These differences are 
reconciled by noting that the radiometer equation assumes 
a uniform antenna density p an t, while using equations pop 
and pip to calculate the rms noise fluctuation in visibil- 
ity space accounts for an arbitrary antenna density. For the 
Pant oc r~ 2 antenna density profile of the MWA, the reduced 
antenna de nsity at large r leads t o a reduced baseline density 
at large U (jBowman et al.ll2006T ). and hence increased noise 
levels at high angular resolution. Figure [3] demonstrates the 
increased noise levels for the largest values of U ~ 800 for 
the MWA operating at v = 165 MHz. In summary, in the 
case of the MWA, using the radiometer equation (and hence 
assuming a uniform antenna density) results in an underes- 
timate of noise for large baselines and an overestimate for 
small baselines. 
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Figure 4. Left: Sample mock spectra of quasar HII regions using the MWA showing the recovered step function (dashed). Center and 
Right: Marginalized probability functions of recovered neutral fraction xjji and radius R for an ensemble of simulations (grey). One 
example (bold) corresponds to the mock spectrum in the left-hand panel. Likelihood functions £(xjji|(^hi}) f° r an observed neutral 
fraction xm are shown by the dashed curves. The vertical lines show the global neutral fraction of the surrounding IGM and the modeled 
region radius R q . Six cases are shown corresponding to tj nt = 100 hr and 1000 hr when (^hi) ~ 0.15, 0.3 and 0.5 (top to bottom). The 
same IGM realization and system noise (scaled by 1/ V^iiit) have been used at each redshift to facilitate comparison between integration 
times. 
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4.2 Continuum foreground subtraction 



to-noise to detect an HII region we fit a step function 



Foreground contamination and its removal can have signif- 
icant consequences for the detectability of the weak, red- 
shifted 21 cm cosmic signal. Both spatial and frequency char- 
acteristics of foreground emission will affect our ability to 
subtract these contaminants sufficiently in degree and accu- 
racy. For the remainder of this paper we concern ourselves 
with the detection of HII regions surrounding known high 
redshift quasars in the spectral regime. We defer discussion 
of the detectability of such HII regions in th e full image 
regim e using a filter matching technique (e.g. iDatta et al.l 
120071 ) to future work. 

Foreground contamination will be brighter than the cos- 
mological 21 cm signal by at least four orders of magni- 
tude. We assume that the total foreground si gnal, which is 
domi nated by galactic synchrotron radiation l|Shaver et al.l 
Il999h . has had extragalactic point sources removed, leav- 
ing a spectrally smooth foreground sign al that will be re- 
movable through contin uum subtraction l|Gnedin fc Shaveil 
120041 ; IWang et alj|2006h . 

The process of continuum foreground removal will 
remove any cosmological signal for line of sight modes 
with fcn <C 2n/y, w here y is the depth of the survey 
l|McQ.uinn et al.ll2006l ). Since HII regions may have sizes cor- 
responding to a significant fraction of the frequency band 
pass, the insensitivity to large scale modes may hamper the 
detectability of HII regions. The MWA will operate with a 
processed bandwidth of 32 MHz. For comparison, the fre- 
quency domain of our synthetic spectra is Av <; 20 MHz. 
Thus, foreground modes along the line of sight with wave- 
lengths longer than the 32 MHz band pass are not present in 
our mock spectra, with the exception of an overall constant 
corresponding to the DC mode of the spectrum. 

In this paper we make the assumption that foregrounds 
could be removed over the 20 MHz of our simulated spec- 
trum, given access to the 32 MHz of spectrum that will be 
available for the MWA. 



4.3 Estimate of signal-to-noise ratio in detection 
of neutral IGM surrounding HII regions 

We next present mock spectra of HII regions surrounding 
known high redshift quasars at times when the global neu- 
tral fraction is (arm) ~ 0.15, 0.3 and 0.5, and discuss the 
prospects for their detection. 

Line of sight spectra are simulated using a synthesized 
beam of width # D = 3.2 '. Figure 2] shows sample mock spec- 
tra for 100 hr and 1000 hr integrations centered on 186, 175 
and 164 MHz (corresponding to z = 6.7, 7.1 and 7.65 re- 
spectively). Instrumental thermal noise in these spectra is 
generated following the prescription described in § 14.11 The 
presence of an HII region could be inferred from this data 
via the detection of a brightness temperature contrast be- 
tween the ionized IGM within the HII region and the neutral 
hydrogen outside. We do not know a priori the radius of the 
HII region. However from observations of the proximity zone 
surrounding a high redshift quasar we know that if an HII 
region is present, its size must be larger than the redshift 
corresponding to the observed transmission blueward of the 
Lya line in near-IR spectra. 

To test whether the mock spectra are of sufficient signal- 



Tmodci(^) = Tjgm — AT x Q(R r , 



R) (13) 



which has three free parameters. These parameters cor- 
respond to the brightness temperature outside the region 
Tigm, a radius i? rcg ion (which we assume to have a value 
R > i? q ) and a temperature brightness contrast across the 
HII region boundary AT. The parameter Tigm includes the 
contribution of the DC component of foregrounds which is 
therefore removed during the fitting procedure. We note that 
our fit is made to the spectrum in only one pixel of the syn- 
thetic map, which includes contributions of foregrounds from 
all regions within the synthesized beam. Using the calculated 
thermal noise ot in each of JV cri channels for the MWA with 
an assumed integration time tint , we then calculate the value 
of x 2 , where 



E 



T% T n io del 



(Tt 



(14) 



for each set of parameters, and hence the likelihood 

£(AT, Region, Tgm) = exp (-x 2 /2) . (15) 

From Bayes' theorem the a posteriori probability for the 
parameters given the observed spectrum follows from 

d 3 P 

dATdft-ogiondTlGM 

OC C(AT, -Rregion, TlGM 



HP ■ HP ■ HP ■ 

Ui1 prior prior prior 
dAT (IRrcgion dTlGM 



(16) 



where we assume prior probability distributions 



dAT 



Ct-Rreg 

dP u 



Q{Rr, 







Ri, 



rfTiGM ^ ^ 

Here is the Heaviside step function, i? m i n = R q and 
-Rmax = 15 proper Mpc. The value of Timax represents an 
upper limit on the size of a proximity zone within which the 
quasar ionizing flux is dominant. If the HII region extends 
beyond i? ma x along a particular line of sight, we effectively 
count this ionization as part of the general two phase IGM. 
Rmin is set by the observed radius of Lya transmission (indi- 
cating a highly ionized IGM) . Provided the quasar is beamed 
with an opening angle > a, where a = arctan(i?b/-R q ), this 
assumption will be valid. This is likely since at z « 6.5 with a 
synthesized beamwidth of ~ 4', a « 13°. We then marginal- 
ize the probability distribution to obtain constraints on the 
parameters of interest (AT and -Rrcgion). Finally, the neu- 
tral fraction is linearly related to the retrieved temperature 
contrast through arm = AT/{22 [(1 + z)/7.5] 1/2 }mK. 

The uncertainty in an observation should arise from a 
combination of IGM fluctuations and instrumental noise. 
This instrumental noise is thermal and therefore Gaussian, 
whereas IG M fl uctuations are highly n on-Gaussian (e.g. 
iLidz et all 2007a; IWvithe fc Moralesll2007l ) . A non- Gaussian 
noise component renders the \ 2 statistic suboptimal. How- 
ever, in the case of the MWA, instrumental noise will domi- 
nate in narrow frequency channels which justifies the use of 
X 2 in estimating the likelihood. 
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4-3.1 Recovered parameters from mock spectra 

In the central and right-hand panels of Figure [4] we show 
marginalized a posteriori probability density functions for 
the recovered neutral fraction xm and radius R ieg ion for 
100 hr and 1000 hr integrations. These parameter estimates 
correspond to the mock spectra in the left-hand panels of 
Figure [4] The vertical lines show the global neutral fraction 
of the surrounding IGM and the modeled region radius R q . 

In addition to the telescope noise, the spectra plotted 
in Figure U show that the recovered radius will be sensitive 
to the large fluctuations in brightness temperat ure due to 
the inhomogeneous ionization state of the IGM (|Lidz et al.l 
2007b). In contrast, the estimates of neutral fraction along 
a particular line of sight average over many fluctuations for 
a range of recovered radii. As a result the neutral fraction 
along a particular line of sight should be measured in an 
unbiased way. On the other hand, different lines of sight will 
show a varying systematic offset from the mean due to the 
inhomogeneous reionization of the IGM. 

To investigate the extent of systematic offset due to 
patchy reionization we have analyzed mock spectra from an 
ensemble of HII regions and lines of sight. Probability distri- 
butions for the best fit parameters R and xm are estimated 
from this ensemble and are shown as grey lines in Figure [4] 
These distributions quantify the level of systematic uncer- 
tainty that is present. In cases where the telescope noise is 
small, the systematic offsets may be larger than the statisti- 
cal uncertainty in the neutral fraction. This is demonstrated 
by the 1000 hr examples in Figure [4] The inhomogeneous 
nature of the IGM, combined with system noise conspire to 
make the determination of R very uncertain despite high 
statistical precision. 

We define the confidence for detection of an HII region 
C by the ratio of the mean recovered neutral fraction to the 
standard deviation of the marginalized a posteriori proba- 
bility density functions for the recovered neutral fraction. 
Cumulative probabilities of detection are calculated using 
an ensemble of mock line of sight spectra. Figure [S] shows 
the cumulative probabilities of detection using single line 
of sight spectra for 100 hr and 1000 hr integrations (dark 
lines). These results indicate that it will be possible to reg- 
ularly detect bubbles at high significance [P(C > 5) > 0.5] 
from spectra of individual HII regions using integrations of 
1000 hr, even at redshifts close to overlap when the neutral 
fraction is as low as m 15%. However at neutral fractions of 
> 30%, measurement will be possible in only 100 hr. This 
is despite the fact that at these late times the typical HII 
region becomes comparable in size to the HII regions gener- 
ated by the quasar (at z < 7 in this model). As stressed in 
§[TJ even when the impact of the quasar is no longer domi- 
nant, the method of detecting the 21 cm signal a posteriori 
along the line of sight to a known quasar will still succeed, 
because neutral gas is known a priori to be absent in a large 
volume surrounding the quasar. 



4.4 Extracting neutral fraction from observations 

The discussion in § 14. 31 describes estimation of the certainty 
with which an HII region can be detected. Once the exis- 
tence of neutral gas is established, the interesting question 
becomes whether the HII region can be used to measure the 
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Figure 5. Cumulative probabilities of detection for 100 hr 
(dashed) and 1000 hr (solid) integrations using single (black) and 
stacked spectra (grey) when the globally averaged neutral fraction 
is (xhi) ~ 0.15, 0.3 and 0.5 (top to bottom). 



global neutral fraction. Below we outline how this might be 
achieved. 

The modeling described could be used to generate a 
likelihood function for an observed neutral fraction xm given 
a trial global mean (xm), 

1 N rIP 

^K^» = F E^ ( • (18) 

i=l (SHl) 

where TV is some large number of model spectra and the 
dPi/dxm are the corresponding probability distributions for 
xm (e.g. those shown as grey lines in Figure[4]| . The resulting 
likelihood functions for the example neutral fractions are 
shown in Figure [4] by the dashed curves. These likelihood 
functions give an idea of the uncertainty on a measurement 
of (xhi) based on an observation of a single HII region. In 
the case of a 1000 hr integration, we find the variance of 
these likelihood functions to be o~ xm « 0.07, 0.09 and 0.1 
for (xhi) = 0.15, 0.3 and 0.5 respectively. The corresponding 
relative errors Cxm/ixm) are 0.5, 0.3 and 0.2. These values 
represent estimates of the achievable uncertainty. 

Given an observed spectrum with probability 
dPohs/dxm for xm along that line of sight, it would 
then be possible to find the a posteriori probability for the 
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globally averaged neutral fraction (xm), 



dP 



d(x m ) 



■£{x m \{xm)) r , (19) 



dxi 



d{x m ) ' 



where dP pr i or /d(a;Hi} is the prior probability for the neutral 
fraction along that line of sight. 

4.5 Recovered parameters from stacked spectra of 
a number of HII regions 

In the previous section we presented results for mock obser- 
vations of single HII regions like those expected to be present 
around luminous high redshift quasars. The observed num- 
ber density of quasar s with luminosities gr eater than L scales 
as N(> L) oc L~ 2 (|Richards et al.ll2006h . Thus while cur- 
rent surveys only detect ~ 1 luminous quasar per MWA 
observing field, in the future, deeper surveys will uncover 
less luminous high redshift quasars at greater number den- 
sities. This increase in density towards lower luminosity will 
make it possible to stack multiple spectra, thereby reducing 
fluctuations in the ionization state of the IGM along partic- 
ular lines of sight which cannot be r emoved through a long 
integration time l|Wvithe et al.|[2005h . In this section we in- 
vestigate whether an increased sensitivity for the detection 
of contrast between neutral and ionized regions could be ob- 
tained by stacking spectra of smaller HII regions around a 
number of less luminous high redshift quasars. We consider 
a value of R q ~ 21 comoving Mpc, which may be more ap- 
propriate for a fainter population of quasars than the 34 
comoving Mpc we assume for the single HII regions studied 
in the previous subsection. 

Mock spectra are simulated as before for 10 HII re- 
gions surrounding host halos of minimum mass 10 12 M Q , 
and then stacked to obtain a single averaged spectrum. The 
stacked spectra are more representative of the average IGM 
surrounding quasars and have decreased levels of instrumen- 
tal noise. We find that averaging over the lines of sight to 
many quasars greatly reduces the systematic component of 
uncertainty in the neutral fraction introduced by variation 
among individual lines of sight. We fit step functions to the 
mock stacks of observed spectra with three free parame- 
ters and estimate the probability distributions of best fit 
parameters from an ensemble of simulations as before. Fig- 
ure [S] shows the cumulative probabilities for detection using 
stacked spectra from 100 hr and 1000 hr integrations (grey 
lines). These results indicate that it will be possible to detect 
bubbles at high significance [P(C > 5) > 0.5] from stacked 
spectra, even at redshifts close to overlap and with modest 
integration times of ~ 100 hr. 



5 DISCUSSION 

We have analyzed the prospects for detection of HII regions 
surrounding luminous high redshift quasars in a percolat- 
ing IGM prior to the end of the epoch of reionization. We 
find that fluctuations in the HI distribution of the surround- 
ing IGM will limit the accuracy (though not the precision) 
with which the global neutral fraction may be recovered ow- 
ing to the finite length of the observed spectra. Thus the 
fluctuating nature of the IGM will lower the signal-to-noise 
ratio achievable on the detection of neutral gas in the IGM 



surrounding high redshift quasars. This degradation of the 
signal-to-noise becomes more serious as reionization pro- 
ceeds since the impact of the quasar on the ionization of 
the surrounding IGM becomes progressively less significant. 
At these late times in the reionization epoch, the fluctua- 
tions in the ionization structure of the IGM that accompany 
the percolation of galaxy HII regions will also significantly 
affect the morphology of HII regions surrounding the high 
redshift quasars, leading to the formation of complicated, 
non-spherical boundaries for isotropic sources. 

Nevertheless, even at late times when the global hydro- 
gen neutral fraction is as low as (xhi) ~ 15% the large vol- 
ume of ionized gas which is known a priori to surround the 
quasar will imprint a detectable signature on the observed 
spectrum that is distinct from a region of typical IGM. We 
find that the presence of neutral gas surrounding a single 
HII region centered on a luminous high redshift quasar will 
be detectable in as little as 100 hours {{xm) J> 30%) using 
the MWA in over 50% of cases at a significance greater than 
5cr. A highly significant detection of neutral gas will be pos- 
sible in 100 hr for a stack of 10 smaller 3 proper Mpc HII 
regions. In both cases an accurate determination of the value 
of neutral fraction will be hampered by lack of knowledge 
regarding the detailed shape of the HII region. We estimate 
that fluctuations in the overlapping IGM between different 
lines of sight will limit the accuracy with which the global 
neutral fraction can be inferred from a single HII region to 
~ 50%, 30% and 20% for (xm) w 0.15, 0.3 and 0.5 respec- 
tively. 

The results of this study imply that the detection 
of neutral IGM surrounding high redshift quasars could 
be achieved utilizing either a single luminous quasar, or 
the stacked observations of a number (~ 10) of less lumi- 
nous quasars. Obtaining a sufficient number of high red- 
shift quasars within a single MWA field to allow for the 
possibility of stacking spectra will require an optical/near- 
IR survey for high redshift quasars in the southern sky 
with flux limits around half a magnitude fainter than ex- 
isting northern sky surveys (|Wvithe et al.l l2005h . Several 
suitable surveys for high redshift in the southern sky are 
currently planned (e.g. SkyMappeiQ& VISTA0). The direct 
observation of neutral gas surrounding the HII regions of 
high redshift quasars discovered in these surveys would be 
complimentary to the statistical detection of an ensemble 
of HII regions via meas urement of the power spectrum of 
21 cm fluctuations (e.g. Tozzi et al. 2000l; iFurlanetto et al.l 



2004al; lLoeb fc Zaldarriagal |2004| ; IWvithe fc Morales! 12007 1 
Iliev et al .1 12006). and would represent a direct detection of 
the end of the reionization era. 
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